\(\text{N}_i\) is the number of emails sent for campaign \(i\),
\(\text{y}_i\) is the number of emails clicked for campaign \(i\),
\(\theta_i\) is the click-through rate (CTR) for campaign \(i\).
We assume a Beta(1,1) prior for each \(\theta_i\), which represents no strong prior belief about the CTR:
\(\theta_i \sim \text{Beta}(1, 1)\)
This is a uniform prior on the interval [0, 1].
2. JAGS Model
Here is the corresponding JAGS code for this model:
CTR_model <-"model { # Loop over the three campaigns for (i in 1:n_campaign) { # Likelihood: binomial distribution for the number of clicks y.i[i] ~ dbin(theta.i[i], N.i[i]) # Prior: Beta(1,1) distribution for the CTR (uniform prior) theta.i[i] ~ dbeta(1, 1) }}"
To implement the JAGS model we need data and paramters to monitor. We can also specify initial values if we wish (JAGS will also do this by default).
Data Input for JAGS:
CTR_data <-list(y.i =c(15, 20, 12), # Number of clicks for campaigns A, B, CN.i =c(100, 120, 110), # Number of emails sent for campaigns A, B, Cn_campaign =3)
Initial Values and Parameters to Monitor:
inits <-function() list(theta.i =runif(3)) # Random initial values for CTRsparams <-c("theta.i") # Monitor the CTRs
Putting all of this together and running JAGS
library(rjags)library(R2jags)CTR_model <-"model { # Loop over the three campaigns for (i in 1:n_campaign) { # Likelihood: binomial distribution for the number of clicks y.i[i] ~ dbin(theta.i[i], N.i[i]) # Prior: Beta(1,1) distribution for the CTR (uniform prior) theta.i[i] ~ dbeta(1, 1) }}"CTR_data <-list(y.i =c(15, 20, 12), # Number of clicks for campaigns A, B, CN.i =c(100, 120, 110), # Number of emails sent for campaigns A, B, Cn_campaign =3)inits <-function() list(theta.i =runif(3)) # Random initial values for CTRsparams <-c("theta.i") # Monitor the CTRsmod <-jags(data = CTR_data, inits = inits, parameters.to.save = params, n.iter =4000,n.burnin =2000, model.file =textConnection(CTR_model))
Compiling model graph
Resolving undeclared variables
Allocating nodes
Graph information:
Observed stochastic nodes: 3
Unobserved stochastic nodes: 3
Total graph size: 11
Initializing model
To check for convergence we can have a quick look at the summary output and trace plots.
library(coda)# turn the model into an mcmc objectmod_mcmc <-as.mcmc(mod)# get trace plot and densityplot(mod_mcmc)
3. Interpretation of Results
After running the MCMC sampling (e.g., using 4,000 iterations and discarding the first 2,000 as burn-in), you will obtain samples from the posterior distributions of \(\theta_i\).
(a) Plot the posterior distributions of the CTRs for the three campaigns using a package like tidybayes in R (or coda or the base plotting functions).
library(tidybayes)library(tidyverse)## get the output in matrix formatmcmc.matrix <- mod$BUGSoutput$sims.matrix## get indexes for thetacampaign_index <-1:3## plotsmcmc.matrix %>%spread_draws(theta.i[campaign_index]) %>%ggplot(aes(y =factor(campaign_index, labels =c("A","B", "C")), x = theta.i)) +stat_halfeye() +ylab("Campaign")
(b) Compare the distributions to assess the relative effectiveness of the campaigns. It appears that the posterior distribution for \(\theta_B\) is shifted to the right compared to \(\theta_A\) and \(\theta_C\), suggesting that Campaign B is likely more effective, i.e., the click through rate is estimated to be higher. Note however, there is an overlap in terms of uncertainty in the estimated click through rate.
4. Additional Analysis
(a) Calculate the posterior probability that Campaign B has a higher CTR than Campaign A and Campaign C:
where \(S\) is the number of posterior samples, and \(I(\cdot)\) is an indicator function. You can do similar calculations for \(P(\theta_B > \theta_C)\).